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CONVERGENCE OF FINITE DIFFERENCE SCHEMES TO THE 
ALEKSANDROV SOLUTION OF THE MONGE-AMPERE 

EQUATION 

GERARD AWANOU AND ROMEO AWI 


Abstract. We present a technique for proving convergence to the Aleksandrov 
solution of the Monge-Ampere equation of a stable and consistent finite difference 
scheme. We also require a notion of discrete convexity with a stability property and 
a local equicontinuity property for bounded sequences. 


1. Introduction 


Given an orthogonal lattice with mesh length h on a convex bounded domain hi C 
with boundary dfl, we are interested in convergent finite difference approximations 
of the problem: find a convex function u G C(h2) such that 

det D 2 u = v in Q 

, 0 4 - 1 ) 

u = g on oil, 


where v is a finite Borel measure and g G C(<9h2) can be extended to a convex 
function g G 0(0). When u G U 2 (0), det D 2 u is the determinant of the Hessian 

matrix D 2 u = ( d ^ x ) • In the general case, the expression det D 2 u denotes 

V *’ V i,j=i,...,d 

the Monge-Ampere measure associated with u. 

Let f \ denote the computational domain and dQh its boundary. Let fh > 0 be a 
family of mesh functions which converge to v as measures. We consider the problem 
with unknown a mesh convex function Uh 

h d M h [u h ] = h d f h in n h ^ 

u h = gondVt h . 


Here Aih[ v h] denotes a stable and consistent discretization of det D 2 v for a smooth 
convex function v. There are several notions of discrete convexity. We require that 
the uniform limit on compact subsets of mesh convex functions is a convex function 
and that a locally bounded sequence of such functions is locally equicontinuous. Of 
course we also require (11.21) to have a solution. A sufficient condition is degenerate 
cllipticity and Lipschitz continuity as defined by Oberman [15], We show that a 
family of solutions Uh of (11.21) converges uniformly on compact subsets to the unique 
Aleksandrov solution of (11.11) . 
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The Monge-Ampere equation (11,11) is a fully nonlinear equation which arises in several 
applications of great importance, e.g. optimal transportation and reflector design. 
Problems in affine geometry motivated the study of the Dirichlet problem. 

The equation det D 2 u = v with v a sum of Dirac masses, and with Dirichlet boundary 
condition was solved by Pogorelov [[18]. For the so-called second boundary condition 
we refer to m Chapter V section 3]. When the measure v is absolutely continuous 
with respect to the Lebesgue measure, with density / > 0 and / G C(O), the con¬ 
vergence of a scheme of the type (11.21) was proved in [9] using the notion of viscosity 
solution. In this paper, we use the notion of Aleksandrov solution, the consistency of 
the discretization (jl . 2j) and approximation by smooth functions to handle the Monge- 
Ampere measure. In [2] the notion of Aleksandrov solution was also used along with 
a different procedure for approximation by smooth functions. 

We note that the method introduced in [6] which is more effective when the measure v 
is a combination of Dirac masses is not consistent. Our requirements for convergence 
which are stability, consistency and solvability of (II. ip . as well as stability under 
uniform convergence on compact subsets of discrete convex mesh functions along 
with local equicontinuity of locally bounded sequences of such functions, are met by 
the finite difference scheme introduced in |9j. However our numerical results indicate 
that a very good initial guess is required for an iterative method for solving the 
nonlinear problem (11.21) if one uses the discretization proposed in [9]. In our numerical 
experiments, the discrete problem (j 1.2 [) is solved with a time marching method which 
has also been used in jTj . The difficulty of capturing singular solutions may be related 
to the choice of the method for solving the nonlinear equation (11.21) . We plan to use 
in a subsequent work semi-deterministic algorithms of the type introduced in [14] 
for hyperbolic Monge-Ampere equations. This is motivated by the observation that 
discretizations of the type introduced in [9] may have multiple solutions. 

If v has density / > 0 and / G 0(0), under our assumptions, Uh converges uniformly 
on compact subsets to the unique viscosity solution of (11.11) if the latter is known 
to have a unique viscosity solution. This follows from the equivalence of the notion 
of viscosity and Aleksandrov solutions, the proof of which we outline. If / > 0 and 
/ G 0(0), a continuous viscosity solution of (11.ip is also an Aleksandrov solution of 
m im Proposition 1.7.1]. The result is also valid for / > 0 and / G 0(0). Indeed, 
consider the problems det D 2 u e = f + e, e > 0. By m Lemma 5.1 ], u e converges 
uniformly on compact subsets to u. One then uses the equivalence of the notion of 
viscosity and Aleksandrov solutions in the non degenerate case m Propositions 1.7.1 
and 1.3.4] and the stability of the notion of Aleksandrov and viscosity solutions under 
uniform convergence on compact subsets HH Lemma 1.2.3] and [5] Theorem 2.3]. 

In this paper we provide the convergence proof of a time marching method for solving 
the nonlinear problem (II.2|) . The main contribution of this paper is the method of 
proof for convergence of finite difference schemes satisfying our assumptions. The 
numerical results indicate that such schemes may not lead to an effective numerical 
algorithm. Our results clarify the nature of efficient discretizations for (II.ip . Another 
consequence of our results is the equivalence of the notions of viscosity and Alek¬ 
sandrov solutions for / > 0 and / G O(D) 0 .0(0). Indeed as we show in section 
14.21 Uh obtained through (1 1.2 p and the discretization proposed in [9] converges to the 





















3 


Aleksandrov solution. It is also known to converge to the unique viscosity solution of 
(11.11) when the latter exists. Hence the result. 

The paper is organized as follows: in the next section we give some notation, recall 
key results on the Aleksandrov solution and finite difference schemes. In section [3] 
we prove the claimed convergence result. We conclude with a proof of convergence of 
the time marching method and numerical experiments. 


2. Preliminaries 

In this section, we recall key results on the Aleksandrov solution of the Monge-Ampere 
equation. We then associate discrete measures to mesh functions. For a smooth 
solution of (II.ip we immediately get a discretization of the Monge-Ampere measure. 
Finally, we introduce finite difference schemes. 


2.1. The Monge-Ampere measure. In this paper, we take the analytic approach 
to the Monge-Ampere measure [T9] . Let K(Ll) denote the cone of convex functions 
on fl and let M(fl) denote the set of Borel measures on 9. For v G C 2 (fl) fl K (9), 
and given a Borel set B , we define a Borel measure AA [w] by 

M[v\(B) = / det D 2 v{x) dx. 

J B 

The topology on M(fl) is induced by the weak convergence of measures 

Definition 2.1. A sequence /j, n G M (9) converges to p G M(fl) if and only if 
p n (B) — > p(B) for any Borel set B C 9 with p(dB) = 0. 


We note that there are several equivalent definitions of weak convergence of measures 
which can be found for example in [7, Theorem 1, section 1.9]. We have 

Proposition 2.2. [Ill] Proposition 3.1] The mapping AA maps C{Vt)-bounded subsets 
ofC 2 (Q) fl K(Q) into bounded subsets of M(fl). Moroever M. has a unique extension 
to a continuous operator on K (9). 


For a convex function v, we will refer to AA\v] as the Monge-Ampere measure asso¬ 
ciated with v. It can be shown that it coincides with the notion of Monge-Ampere 
measure obtained through the normal mapping [19, Proposition 3.4]. 


Definition 2.3. Given a Borel measure v on 9, a convex function u G (7(9) is an 
Aleksandrov solution of 

det D 2 u = v, 


if and only if AA[u] = u. 


We recall an existence and uniqueness result for the solution of (11.11) . 

Proposition 2.4 ([T3] Theorem 1.1). Let Q be a bounded convex domain of 
Assume v is a finite Borel measure and g G (7(99) can be extended to a function 
g G (7(9) which is convex in 9. Then the Monge-Ampere equation (11.11) has a unique 
Aleksandrov solution in K(Q) fl C(9). 
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Throughout this paper, we will follow the convention of denoting by p a measure v 
which is absolutely continuous with respect to the Lebesgue measure and with density 

p. 

2.2. Discrete measures associated with mesh functions. Let h be a small po¬ 
sitive parameter and let 

Z^ = { mh, m G Z d }, 

denote the regular uniform grid of M d . By a mesh function we mean a real-valued 
function defined on Z^J. We denote by Ch the cone of discrete convex mesh functions. 
We consider in section 14.21 a notion of discrete convex mesh functions which fulfill the 
assumptions of this paper. 

The computational domain is defined as = hi fl Z^ and its boundary is simply 
dVt h = { x G hi fi Z^, x <£ Qh }■ 

Let Vh be a mesh function such that Vh > 0 on 12^. We associate to Vh a Borel measure 
which we denote here by ry and defined by 

v h (B) = h d v h (x), (2.1) 

xen h nB 


for any Borel set B. 

Given a continuous function v on f2, we use the notation rh(v) to denote its restriction 
to 12/,,. If v > 0 is a continuous function on hi, we have 


lim r h (v){B) = / v(x)dx. 
h ^o J B 

for any B C 12 satisfiying \dB\ =0. In other words, the measures ty converge weakly 
to v. 


Definition 2.5. Let Vh G Ch for each h > 0. We say that Vh converges to a convex 
function v G C(f2) uniformly on compact subsets of LI if and only if for each compact 
set K C 12, each sequence hk —> 0 and for all e > 0, there exists h-± > 0 such that 
for all hk, 0 < hk < h-i, we have 

max \v hk (x) - v(x)\ < e. 
xeKnzf 


2.3. Finite difference schemes. Let > 0 be a family of mesh functions which 
converge to v as measures and let Mlh [ry,u] denotes a discretization of det D 2 v for a 
smooth convex function v. We are interested in the discrete Monge-Ampere equation 

o. 

The discretization Miy,, [rh v \ is said to be consistent if for all C 2 convex functions v, 
and a sequence Xh G Ll h such that Xh —> x G 12, lim^o Mlh[ r h v \{ x h) — det D 2 v(x). 

The discretization is said to be stable if the problem Mih[vh] = fh has a solution Vh 
which is bounded independently of h. 
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3. Convergence 


We recall that our assumptions are that the nonlinear equation (jl.2jl is solvable, the 
discretization is stable and consistent and finally the uniform limit on compact subsets 
of mesh convex functions is a convex function and that a locally bounded sequence 
of such functions is locally equicontinuous. 

Lemma 3.1. There exists a constant Co such that for all Borel sets B C Tl 

\Mh[w h ](B) - M h [vh](B)\ < C 0 max \w h (x) - v h {x)\. 

ze-Bnz^ 


Proof. We have by 02.11) 


\M h [w h ](B) - M h [v h ](B)\ 


h d 


Wk ( x ) 

xGBnn h 


Vh ( x ) 

xGBntth 


< ^ h d \w h (x) - v h {x)\ 

x&BCiClh 


< 


max \w h (x) - v h (x)\ h d 

xeBnzt “ 

xeBnrih 


<C max \w h (x) - v h (x)\. 
xeBn'Z'l 


□ 


The constant Co in Lemma [3. II satisfies 

5^ h d < C 0 . (3.1) 

x 

Lemma 3.2. For v G K(fl) flC 2 (h2), h d M. h [rhv\ converge weakly to det D 2 v. 


Proof. By consistency, M.h\rhV]{x) —>• det D 2 v(x) for v G it (Q) nC 2 (h2). By dehni- 
tion of the Monge-Ampere measure and of discretization of the integral, and using 
consistency, we have for a Borel set B with Ai[v](dB) = 0 


det D 2 v{x) dx 


' B 


Moreover using (13.1|) 


lim h d 

h —^0 


£ 

xGBnflh 


det D 2 v(x). 


h d ^ det D 2 v(x) 

x&Bnci h 


li d M h [r h v}(B ) 


= ^ h d (det D 2 v(x) — Aih[rhV]{x )) 

x^Bnfih 

< Cn max I det D 2 v(x) — Aih\rhv](x)\, 
xcBnn h 


from which the result follows by consistency. 


□ 


We have the following weak convergence result for discrete Monge-Ampere measures 

Theorem 3.3. Let Vh G Ch converge uniformly on compact subsets to v G K(Q) fl 
C(f2). Then h d M.h{vh\ converges weakly to M.[v\. 
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Proof. Let v t G iL(O) nC 2 (0) converge uniformly to v on O. The existence of may 
be proven as in [3|j. Let B be a Borel set with Ml[v](dB) = 0. Given <5 > 0, we seek 
ho > 0 such that \h d Mlh[vh\{B) — Ml[v\(B)\ < 5 for all 0 < h < h 0 . 

By Proposition 12.21 3e 0 > 0 such that \Ml[v eo ](B) — M[v\(B)\ <5/3. By Lemma [3721 
3ho > 0 such that for all 0 < h < ho, \ h d Mlh[(rh(v eo )](B) — Ml[v eo }(B) \ < 5/3. 

We may assume that max^^ |u eo (a;) — u(x)| < 5/(6Co) and since Vh converges to v 
on B, we may assume that for h < h 0 , max xg -g nZ d |n(a;) — Vh{x)\ < 5/(6Co). Thus 
we have max xg ^ nZ d |v eo (x) - v h (x)\ < 5/(3 C 0 ). By Lemma EHJ \M h [r h (v eo )\(B) - 
M.h[vh](B)\ <5/3. This concludes the proof. 

□ 


We can now prove the main result of this paper 

Theorem 3.4. The mesh function Uh defined by ( 11 . 21 ) converges uniformly on compact 
subsets to the Aleksandrov solution u of (11.11) . 

Proof. By the stability assumption, the family Uh is uniformly bounded and by our 
assumption on the discretization, locally equicontinuous. By the Arzela-Ascoli theo¬ 
rem, there exists a subsequence Uh k which converges uniformly on compact subsets to 
a function v. Since Uh £ Ch the function v is convex by our assumptions on discrete 
convex functions. Since Uh is uniformly bounded, v is convex and bounded on fl, 
hence continuous on Arguing as in the proof of [3( Theorem 4.3] one proves that 
v G C(f2). Since h d A4h[ v h] converge weakly to Ml [u] and Uh = g on dQ, the function 
v is an Aleksandrov solution of (11.11) . By uniqueness, v = u and hence the whole 
family Uh converges uniformly on compact subsets to u. □ 


4. Convergence of a time marching iterative method 

Let us denote by Mt(Q h ) the set of mesh functions, i.e. the set of real valued functions 
defined on Q h . Since flh is a finite set, there is a canonical identification of A4(h2 /l ) 
with for some integer N. We will now also use the restriction operator for 
vector and matrix fields. For x G M. N , |x| = (^2f =1 xj)2 denotes the Euclidean norm 
of x and |x|oo = maxj=i ) ... i jv |xj| denotes its maximum norm. 

We make the assumption that the mapping Mlh is Lipschitz continuous with Lipschitz 
constant K > 0 i.e. 


I M h [v h ] - MlhlwhWoo < K\v h - Whlo^v^Wh G R n . 

Here we make the abuse of notation of identifying a mesh function with its vector 
representation. We also make the assumption that problem (1 1.2 j) has a unique solution 
Uh which can be computed by a time marching method 

u k h +l = u k h + -M h [u k h \ in 
h 

Uh = r h {g) on (90ft, 


(4.1) 












7 


for n > /x 0 and a suitable initial guess. Such assumptions are satisfied by proper 
Lipschitz continuous degenerate elliptic schemes as defined by Oberman [IS] . Al¬ 
though our theory indicates convergence of the discretization for the case where the 
measure fi is a combination of Dirac masses, we were not able to get numerical evi¬ 
dence of convergence for the above iterative method for the discretization proposed 
in p] even if we use the exact solution as initial guess. Similar results for Newton’s 
method were reported in [6]. 

Let us denote by A h the standard finite difference discretization of the Laplace ope¬ 
rator and let e t denote the i th vector of the canonical basis of M d . For x G 14 and 
v h G A4(l4), we have 

v h (x + hei) - 2v h (x) + v h (x - hei) 

A h v h (x ) = 2^-T2-• 

i— 1 

When the measure // is a combination of Dirac masses we obtained better numerical 
results with the preconditioned iterative method 

- A,.u ‘ +1 = -A*u* + -M h \4} in Q h 
Uh = r h {g ) on <914, 

for n > /i] under the above assumptions. Moreover numerical experiments indicate 
that the method (14.21) converges faster than (14.11) . The idea to use the Laplacian for 
faster iterative methods has a long story in various contexts pj p. 58, and a remark in 
that direction for proper Lipschitz continuous degenerate elliptic schemes was made 
in ra- See also [IB]. We use the terminology preconditioned iterative method for 
(14.2p by analogy with preconditioned techniques for linear equations. An advantage 
of the preconditioned iterative method (14. 2 j) is that fast Poisson solvers and standard 
multigrid methods can be used at each step. 

The proof of convergence of the iterative method (14.21) does not follow the approach 
in [2] for proving convergence of the basic iterative method (14. ip . The proof of the 
latter does not seem to extend to the preconditioned version (14.2p . We take a different 
approach which consists in using the fact that (14.11) converges to the discrete solution 
of (ll.2p and properties of the inverse of the operator A h . 


4.1. Convergence of the preconditioned iterative method. It can be shown 
[12j Theorem 4.4.1, that for / G C(12) the problem 

A h [z h ] = r h (f) infi h 

Zh = 0 on <914, 


has a unique solution. We denote by A h 1 the inverse of the operator Ah with homo¬ 
geneous boundary conditions. Let 1111 denote the operator norm of A^ 1 , i.e. 


A ^ 1 


sup 

\ v h looT^O 


\^h |oo 


By Theorem 4.4.1 of [12], HA^H is bounded independently of h. We note that 
Theorem 4.4.1 of [T2] is proven for dimension n = 2 but the proof extends immediately 
to arbitrary dimension. 
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The main result of this section is the following theorem 

Theorem 4.1. Let Aih denote a Lipschitz continuous finite difference scheme such 
that the mapping Tf : Miflfi) —> M(Ll h ) defined by 

Tfivh] =v h + -M h [v h \, 

is a strict contraction for p > p 0 > 0. Then for some p± > 0, the mapping T 2 : 
A4(O h ) —> M(Ll h ) defined by 

T 2 [v h ] =v h - -A^Mhivh], 

p 

is also a strict contraction for p> p\. 

Proof. By assumption, there exists a constant C i such that 0 < C\ < 1 and 

| l 1 \Ph\ ^ 1 [^h\ | oo A C11 'Wh |ooj 

for all Vh,Wh G M.{Pt h ). One may decompose T 2 — T 2 [wf\ as 

T 2 [v h \ - T 2 [w h ] = T 2 [v h ] - Tfivh] + Tfivh] - Tfiwh] + Tfiwh] - T 2 [w h } 

= (Tfivh] - TiK]) + (T 2 [v h ] - Tfivh]) - (' T 2 [w h \ - !>,]). 

Moreover 

Tfivh] - T 2 [v h ] = - (M h [v h \ + A ~f 1 Mh[vh\) = - (I + A^ 1 ) M h [v h \, 

p p 

where I denotes the identity operator on Ai(Ll h ). We then get 

{T 2 [v h ] - T\[v h \) - (T 2 [w h ] - Ti[w h }) = *-(/ + A^ 1 ) (M h [v h ] - M h [w h ]). 

T 

We recall that M.h is Lipchitz continuous, i.e. 

I M h [v h ] - M h [w h \ |oo < K\v h - iu/jIoq, Vv h ,w h G 
One deduces that 

I T 2 [v h \ - TfiwhWoo < \Tfivh] - T^WhWoo + |- (/ + A” 1 ) (M h [v h ] - Mh[wh])\oo 

T 

K 

< CfiVh - Whloc + —II I + A" 1 !! \v h - w h loo 
A + — \I + A d 1 ||^ | Vh — Wh\oo- 

Since ||/ + A^" 1 1| < ||/|| + || A^ 1 1| is bounded independently of the discretization step 
h and 0 < C\ < 1, one may choose p big enough such that 

CiT-HJ + A- 1 !! < 1, 

p 

making T 2 a strict contraction mapping. This concludes the proof. □ 


Under the assumption of the above theorem, both the iterative methods (14.111 and 
(14.2ft converge linearly to the unique solution Uh of (11.211 . 
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h 


b 

1/2 3 

1/2 4 

1/2 5 

1/2 6 

1/2 7 

1 /2 s 

50 

4.71 10" 1 

2.86 10” 1 

1.69 10" 1 

9.77 lO" 2 

5.50 10~ 2 

3.02 10~ 2 


Table 1. 


4.2. A numerical example. We say that a mesh function Vh is discrete convex if 
and only if A e Vh(x) = Vh{x + e) —2 Vh(x) + Vh(x — e) > 0 for all x G Vth and e G Zf for 
which A e Vh{x) is defined. Then the uniform limit of discrete convex mesh functions 
is convex P, Lemma 2.11]. Moreover a bounded sequence of such functions is locally 
equicontinuous pj. 

Following [9], we define 


M h [v h ](x) 


inf 

(“1 ,-,an)£W h (x) 



V h (x + ai) 


2v h (x) +v h (x-a>i) \ 

N 2 ’ )' 


where for x E flh, Wh(x) denotes the set of orthogonal bases of such that for 
(aq, ..., a n ) E Wh(x), x ± cq E f4, for all i. 

It is known that M.h[vh\ satisfies the assumptions of degenerate ellipticity and Lip- 
schitz continuity as defined by Oberman [T5j. The consistency of the scheme was 
proved in [9] while for a proof of stability, we refer to [3j. Note that M.h[ v h\ A 0 im¬ 
plies that Vh is discrete convex. Hence the discrete convexity assumption is enforced in 
the discretization. Moreover, as pointed out in |3], if one considers Mh[vh](x) +eVh(x) 
where e is taken close to machine precision, the discretization is proper and hence 
uniqueness holds. 


For the numerical experiments, the space dimension d is taken as 2 and the computa¬ 
tional domain is the unit square (0, l) 2 . Numerical experiments with v a Dirac mass 
was reported in earlier papers, e.g. [9]. Here we consider the example of [6] where v 
is the sum of two Dirac masses, i.e. we take 


u(x,y) = 


I y- 


nun 


0 - i) 2 + iv - 1) 2 > J(x - 1) 2 + (y - 1) 2 


if -j < x < 
otherwise, 


and v = 7r/2 <5o/ 4il / 2 ) + 7 t/2 A( 3 / 4 ,i/ 2 ) - For simplicity, we only use a 17 point stencil. 
The initial guess is taken as the exact solution and the nonlinear equations solved 
with (14.20 . Errors are given in the maximum norm and reported on Table [1] 
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